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ABSTRACT 



We present model simulations of time-variable infrared (IR) emission from dust as a consequence of variability of the incident 
radiation. For that we introduce a generalized treatment for temperature variations in a dusty environment, which is not limited to any 
specific astronomical source. The treatment has been incorporated into a simplified clumpy torus model, with the radial brightness 
distribution as the main parameter, to study the IR emission of type 1 active galactic nuclei (AGN). We show that any variability 
signal in the optical is smoothened stronger if the brightness distribution is very extended, and this smoothing strongly depends on 
wavelength. This also affects time lags between the optical and near-/mid-IR emission, which can be up to 10s of sublimation radii 
for long wavelengths and extended brightness distributions. The dependence of time lag on wavelength and distribution can be used 
to quantify the brightness distribution in an AGN torus, either by comparing optical light curves to near-IR and mid-IR light curves, 
or by directly comparing near-IR to mid-IR light curves. Moreover, our model has been applied to near-IR data of the nearby Seyfert 
1 galaxy NGC 4151. We show that the simple model can reproduce the overall observed variability signal and found that about 40% 
of the energy in the variability signal in the l/-band has been converted into AT-band variability. This low value may be explained by 
a "snowball" model of gradually-sublimating clouds at the inner edge of the torus. We also note that our modeling does not support a 
change of time lag/sublimation radius over the observed light curve epoch in spite of a significant change in V-band emission. 

Key words. Galaxies: Seyfert - Galaxies: nuclei - Galaxies:active - Radiative transfer - Galaxies: individual: NGC4151 



1. Introduction 

An optically- and geometrically-thick, circumnuclear dusty re- 
gion (= "dust torus") is one of the key ingre dients of the unifi ca- 
tion scheme of active galactic nuclei (AGN) dAntonuccil 1993b . It 
is held responsible for angle-dependent obscuration of the cen- 
tral accretion disk and broad-line region, explaining the apparent 
difference between type 1 (unobscured line-of-sight) and type 2 
(obscured line-of-sight) objects. Its dust content absorbs part of 
the accretion disk's UV/optical emission and reemits in the in- 
frared (IR). Since the torus resides on parsec scales, direct ob- 
servations have been challenging since its angular size is much 
smaller than the resolution power of single telescopes in the IR. 

With the advent of IR long-baseline interferometry and 
reverberation mapping, it is now possible to det ermine the 
sizes of the dusty region around nearby AGN (e.g. iJaffe et alJ 



and sizes are in general agreement with theoretical expecta- 
tion of the sublimation radius based on thermal equilibrium 
of l arge graphite grains or anisotropic accretion disk emission 
(e.g. IKishimoto et alJ2007l) . These results were confirmed by di- 
rectly mea suring near-IR /T-band ring radii using the Keck inter- 
ferometer ( Kishim oto et al.l l2009bl 1201 lab . It should be noted 
that variability of the accretion disk occurs on multiple time 
scales (from intraday variability up to many years) and ampli- 
tude ranges. The torus has been observed to react on va riabil- 
ity in the range of few days to as long as decades (e.g.j Glass 
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2004 



Tristra m et alll2007t B eckett et al. 2008; Kishimo toet alJ 



2009b; Bu rtscher et alj |2009; Pott et al.ll2010 ). Of special inter- 



est are observations of type 1 AGN where the torus is presum- 
ably seen closer to face-on than in type 2 A GN without too 
much complic ations due obscuration effects dKishimoto et alJ 
2009al I2011hh . While sizes are the primary parameter to ex- 
tract from interferometry data, reverberation mapping uses the 
fact that any variability in the accretion disk emission causes 
a delayed response in the near-infrared continuum emission 
from the hottest, inner region of the dusty torus because of 
the light-travel time from the disk to the torus. Both inter- 
ferometry and reverberation mapping found that the time lags 



IR interferometry has proven as a tool to access more de- 
tailed information about the dust torus than just sizes. B ased 
on mid-IR interferometric observations Jaffe et a D (12004 ar- 
gued that the dust in the torus is confined to small clumps 
rather than smoothly distributed, as theo retically and obser- 
vationally suggested in ear ly studies (e.g. iKrolik & Begelman 
1988; iTacconi et al.l Il994l) . In the following clumpy torus 
models have been successful in reproducing spectral (e.g. 
Nenkova et all 120021; iPolletta et all 120081: iRamos Almeida etatl 
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| 201 It lAlonso Herrero et al.11201 lh and spati al infor mation (e.g. 
Honig et al.l l2006i 120081: ISchartmann et al.l 120081: iHonig et al.l 
2010). Moreover, we recently showed that the dependency of 
visibility on wavelength and spatial frequency can be used 
to constrain the radia l brightness distr i bution of the dust 
(IKishimoto et alJl2009bt IHonig et al.ll20lol:lHonig & Kishimoto! 
2010; IKishimoto et al.l 1201 lbl) . Using such an analysis, it has 
been found that the brightness distribution and the slope of 
the mid-IR spectra are strongly correlated. Both parameters are 
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probably tracers for the radial distribution of the dust, so that 
combining interferometry and photometry provides a tool to 
characterize the properties of the dust torus in different objects. 

At present IR interferometry has some tight limitations re- 
garding brightness and spatial resolution, so that only the bright- 
est nearby sources can be studies in detail. Reverberation map- 
ping, on the other hand, may serve as a complimentary tool for 
objects that are not accessible by interferometry. In this paper 
we present a correspondence in reverberation mapping to the 
interferometric method of determining the brightness distribu- 
tion in type 1 AGN. For that we use a simplified version of a 
clumpy torus mode l, which has been successfully applied to in- 
terferometric data (iKishimoto et al.l2009bl.l201 lbl) . In Sect.|2]we 
describe a theoretical approach to handle temperature variations 
in a d usty medium, generalizing the seminal work of Barvainis 
(119921) . Based upon this theory, we outline a simple model for 
type 1 AGN in Sect.[3]and present model light curves and cross 
correlation function at IR wavelength that show the dependence 
on the brightness distribution of the torus. As a proof of con- 
cept we use our model to reproduce the Zf-band light curve of 
NGC 4151 in Sect. [4] The results are summarized in Sect. [5] 

2. Luminosity and temperature variations in a dusty 
medium 

In this section we describe theoretically how variability in the lu- 
minosity of the incident radiation changes the temperature of the 
dust. The new dust temperatures can then be used to calculate the 
changed emission of the dust. Under the assumption that heating 
and cooling times are negligible, no iterations are required. 

For dust grains in radiative equilibrium of size a at distance 
r from a radiating source with spectrum L v , the received power 
per dust grain, L rec = J L v /(4nr 2 ) g abs; v na 2 dv = L abs a 2 /(4r 2 ) 
equals the emitted power L em = J 4-na 2 (2abs;v nB v {T)dv, where 
2abs;v is the absorption efficiency of the dust and B V {T) denotes 
the Planck function of a black-body at temperature T. Equating 
L rec and L em leads to 

L, bs = I6nr 2 Q ahs ?(T) cr SB T 4 , (1) 
where J 2abs;v xB y (T) dv has been replaced by the Planck mean 

absorption efficiency Qabs;p(T) • <tsbT 4 . If the emission source 
is varying with time, i.e. dL/dt + 0, the temperature of the dust 
grain will change. From Eq.[T]we find 



d^abs 4L a bs 



abs 



abs;P 



(2) 



dT T £ abs;P 

Here, we define Q' bs:P = dQ d \, s -p/dT as the ^-derivative of the 
Planck mean absorption efficiency of the dust. Both the Planck 
mean absorption efficiency and its derivative can be calculated 
from the optical properties of dust and only depend on the dust 
temperature. 

Absorption efficiencies of astronomical dust species 
(graphite, silicates) and compositions typically peak in the 
UV/optical in about the same v-range where many astronomical 
sources radiate (e.g. stars or AGN accretion disks). Therefore 
we assume for practical purpose^] dL a b s /L abs ~ dL/L (with 
L — J L v dv) and obtain the variability function 

, T dL /4 eynr 1 

1 Without this assumption, we would have to replace dL^/L^ = 
dL/L ■ gab S; v/2ab S ;L where j2abs;L = / £v2abs;vdv / / L v dv is the absorp- 
tion efficiency averaged over the emission profile of the source. 
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Fig. 1: Dependence of the relative temperature change dT/T on 
the actual dust temperature for a variability in luminosity of 
dL/L - 1 . The black-solid line shows standard ISM (see text for 
details) while the blue-dashed line is for large graphite grains. 
The red-dotted line indicates the black-body limit. 



For a given temperature T of the dust (e.g. as obtained by 
equilibrium calculations) the change of temperature dT can 
be analytically calculated for any change in incident radiation. 
However the change in temperature depends on the temperature 
itself. On the other hand this additional T-dependence is mi- 
nor. In Fig. Q] we show the dependence of the fractional change 
of temperature dT/T on the dust temperature T for a constant 
dL/L = 1 . The black-solid line shows this dependence for stan- 
dard ISM dust with 53% silic ates, 47% graphite grains (optical 
constants from Draine 120031) and a grain size distribution ac- 
cording to Mathis et alJ ( 19771 MRN). For temperatures below 
~100K and above 1 000 K the fraction is nearly constant at about 
0.16-0.19. Between 100 K and 1000 K the Planck mean opacity 
has a "knee" that is reflected in dT/T as a bump up to 0.24 at 
400 K. The blue-dashed line in Fig.[T]shows large graphite grains 
with sizes 0.07 /vm < a < 1 /mi (MRN size distribution). Below 
100 K the graphite grain dT/T wiggles around the same mean 
value of 0.17 as the standard ISM dust. From 100 K to the sub- 
limation temperature (~1500-2000K) it increases to about 0.2. 
Overall these changes are not very big and may be approximated 
by a constant factor. In fact, if we assume the black-body limit 



abs;P 



0) we get 



dT 
T 



1 dL 
4 T 



(4) 



where the relative change in temperature is directly linked to the 
relative change of luminosity without any extra T-dependence, 
and the constant factor is 0.25. This is shown as a red-dashed 
line in Fig.Q] 

It should be emphasized that this treatment is not restricted 
to a special class of objects (e.g. AGN, see Sect. [3) but can be 
used for dusty media in general. The presented treatment of tem- 
perature variations implicitly assumes that the reemission occurs 
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in the optically-thin wavelength regime. Otherwise an additional 
"self-heating term" has to be included on the left side of Eq. Q] 
and the treatment becomes iterative. However, even in case of 
a clumpy dust torus, where the primary emitting regions in the 
near- and mid-IR are the directly illuminated, optically-thin lay- 
ers of otherwise optic ally- thick dust clouds, the temperature and 
emissio n profile is dominated by d irect heating from the central 
source dHonig & Kishimotol 12010) and the presented treatment 
should be approximately applicable. 



3. IR variability in AGN 

3. 1 . Type 1 variability model 

One of the applications of a variability model can be reverber- 
ation mapping measurements of AGN. Fluctuations in the inci- 
dent radiation from the accretion disk, which primarily emits at 
UV/optical wavelengths, causes th e re-emissio n of the dust in 
the torus to change over time (e.g. Glass 2004). Since the dust 
is located at some distance from the accretion disk, the variabil- 
ity in the IR is delayed with respect to the optical emission by 
the light-travel time from the disk to the dust. This method has 
been exploited to measure the near-IR radius in nearby AGN 
(e.g. lGlassll2004tlSuganuma et al.ll2006l) . 

Reverberation mapping is usually done in type 1 AGN where 
obscuration effects are, in general, significantly lower than in 
type 2 AGN. In these objects the dust-/brightness distribution is 
projected onto a "disk" where the hottest dust e mits at the small- 
est ra dii and the cooler dust at larger distances. iKishimoto et al.l 
(2009a) presented a simple analytic model for the IR emission, 
which has been shown to be a good representation of a clumpy 
dust torus in type 1 AGN dHonig & K ishimoto l2010h . In this 
concept the dust emission S v (r) (dominated by the optically-thin 
layer of an otherwise optically-thick dust cloud; see Sect. [2]) at a 
distance r from the AGN is characterized by the corresponding 
radiative equilibrium temperature T. The total emission L(tor) v 
of the torus is then calculated by inte grating from the sublima- 
tion radius to a (well selected; see Honig & Kishimoto 2010) 
outer radius r out , and accounting for a radial dust distribution 
that is parametrized as a power law with index a: 

L(tor) v = 4tt f M nB v {T{r)) (— ) dr (5) 

In this concept the term (r/r su b) a is directly related to the sur- 
face filling factor of the torus (for details see the discussion in 
iHQnig & Kishimotoll2010l Sect. 2.2). The dust equilibrium tem- 
perature can be obtained either by approximating a power-law 
(iBarvainislI 19871) or self-consistently from Eq. Q] (and using the 
sublimation radius and temperature as reference values) by solv- 
ing the following equation: 

6abs ; pCO / J_\ 4 = /^\ 2 
2abs;p(r su b) \ * sub / \ f f 

For practical purposes this equation can be solved using a pre- 
calculated look-up table for the left term and interpolating for 
the respective r. 

With these prerequisites it is possible to calculate 
wavelength-dependent light curves of type 1 AGN based on 
Eq. |3] considering the traveling time of any variability through 
the torus. 
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Fig. 2: Model light curves at 2.2 yum (red) and 8.5 //m (blue) for 
different values of the radial brightness distribution power law 
index a. From dark to light the distribution changes from very 
extended to very compact. The black-dashed line illustrates the 
variability signal from the accretion disk. 



3.2. Model results 

The model presented in Sec. l3.1l has only one parameter, a, that 
represents the radial brightness distribution of the torus. We have 
shown that this parameter is well correlated with the actual radial 
distribution of th e dust at least in type 1 AG N even in more com- 
plicated models (Honig & Kishimoto 2010). In consequence the 
radial dust distribution can be probed by observables such as the 
mid-IR spectral index or th e comparison between mid-IR and 
near-IR interferometric size (Honi g et al.1 2010; Kishimoto et al. 
1201 lbh . Here we analyze the effect of radial brightness distribu- 
tion (and, in turn, radial dust distribution) on IR light curves of 
AGN. 

3.2.1. Light curves 

In Fig. [2] we show light curves A/(f)//o (where /o = f(to) 
is the flux at time fo, the starting time of the monitoring, and 
A/(0 = /(0 _ /( f o) is the flux difference between to and t) at 
2.2 /im and 8.5 //m for various radial brightness distributions, 
parametrized by the power-law index a. The initial AGN vari- 
ability signal is a step function with a width of 0.5r su b/c and 
an amplitude of 0.5 AL/L, and the light curves essentially show 
the transfer function of this signal. The lighter colors correspond 
to more compact distributions. In these cases most of the light 
is coming from a region close to the sublimation radius, i.e. the 
dust is concentrated to the inner torus. The darker colors have 
more extended brightness distributions. Here a significant frac- 
tion of the total dust mass is located at larger radii and directly 
heated by the AGN. 

We plot the light curves as a function of the intrinsic time 
scale f su b = r m b/c, which corresponds to the light-travel time 
from the source to the sublimation radius r su b. In this way the 
variability signal from the AGN can be considered a tomo- 
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Fig. 3: Model cross correlation function (CCF) for light curves 
of AGN variability at 2.2 /mi (red) and 8.5 /mi (blue). From dark 
to light the distribution changes from very extended (a = 0.0) to 
very compact (a = -2.0). The dashed line marks a lag time of 
r S ub/c 



graphic device to map out the brightness distribution in the torus, 
and elapsed time and distance from the AGN become equivalent 
independent of the object's luminosity (f su b = r su b/c oc L 1 ^ 2 , see 
Sect.EJ. 

The near-IR (2.2 yum) light curves in Fig. [2] show that the 
more concentrated distributions are stronger peaked, while the 
distributions with a closer to have longer tails. In the very 
compact objects almost all of the near-IR light comes from the 
peak black-bodjQ emission of the hottest dust. When the dust 
distribution becomes more extended, then the Wien tails of the 
cooler dust emission start to contribute relatively more since the 
relative amount of hot dust decreases. 

The behavior in the near-IR has its correspondence at longer 
wavelength. The mid-IR (8.5 /mi) light curves with a compact 
brightness profile are peaked in the inner torus where the hot dust 
is located. Over time the brightness decays quickly, however not 
as fast as in the near-IR. As the brightness distribution becomes 
shallower (more extended), the initial peak at r su b/c vanishes, the 
light curves get flatter, and for a > -0.5 they keep on increasing 
out to > 10 x r su b/c. The difference between the mid-IR and the 
near-IR is that in the mid-IR the emission comes predominantly 
from either the Rayleigh- Jeans tail of hot dust emission in the 
inner torus or from the peak black-body emission of cooler dust 
at larger distances. In the near-IR the flux originates predomi- 
nantly in the black-body peak of hot dust with little contribution 
from the steeply dropping Wien tail of cooler dust. For steep 
brightness profiles, the hot-dust Rayleigh- Jeans tail dominates 
the mid-IR emission due to the lack of extended dust. Around 
a « — 1, the contribution of the black-body peak from extended 

2 We note that the correct description of the dust emission is a gray- 
body emission with source function S r = S,,(l - exp(-r v )). For the 
purpose of qualitatively describing the results, however, we will use the 
term "black body". 



dust becomes about equal to the hot dust Rayleigh- Jeans tail, 
which takes over for even shallower brightness profiles. 

3.2.2. Lag times 

Since the contribution from different emission regions changes 
as a function of wavelength and brightness distribution, we can 
expect delays between the peak of the emission in light-curves. 
Classical reverberation mapping tries to quantify the lag time 
between the AGN variation and its correspondence in the near- 
IR by means of the cross correlation function (CCF) 

CCF = J f x (t)f Y (t-At)dt (7) 

where fx(t) and /y(f) are light curves at wavebands X and Y, 
and At is the introduced time lag between both light curves. This 
method can be applied, in principle, also to the mid-IR. We note 
that in the framework of this paper we define lag times as the 
peak in the CCF. 

In Fig.[3]we show the CCF of the AGN variability signal with 
the 2.2 /mi and 8.5 //m light curves, respectively. As expected 
the time lag between AGN variability and near-IR emission is 
close to r SU b/c and the CCF is quite narrow. The fact that it is 
actually slightly larger than unity is in part a result of the used 
variability function (step function with width l/2f su b) and the 
smoothened-out response by the dust. It is seen that the lag time 
depends slightly on a, so that extended emission profiles show 
slightly longer lag times (~ 1 .5 x r su b/c) while the most compact 
brightness distributions are at ~ 1 . 1 x r SU b /c. 

The mid-IR light curves show a different behavior. If the 
brightness distribution is compact, the peak in the CCF is still 
well-defined and close to r su b/c. However when the distribu- 
tion becomes more extended, the CCF becomes very broad. For 
a < -0.5 the peak is no longer well defined and shifts to very 
long time lags. This is expected from the light curves and re- 
flects the change in contribution from hot and cooler dust to the 
mid-IR emission (see Sec. 13.2.1] ). 

The change in CCF (or light curve) with a may be used 
to distinguish between compact and extended dust distributions. 
Interestingly, since the near-IR CCFs and light curves are peaked 
around r su \,/c for all a, it is not necessary to take the AGN vari- 
ability signal as a reference, but instead use the near-IR as the 
CCF-reference signal for the mid-IR. While the actual size of 
r su b cannot be determined in that way, the dust distribution rel- 
ative to r su b is still accessible. In Fig. [4] we show the CCF be- 
tween near-IR and mid-IR light curves. For each a the near-IR 
light curve was used as reference for the corresponding mid- 
IR light curve. As expected the compact brightness distributions 
have no time lag between near- and mid-IR, illustrating the fact 
that the mid-IR emission originates from the Rayleigh-Jeans tail 
of the hot dust emission. With a < -1 the CCF curves become 
significantly wider and the lag shifts to longer time lags. For 
0.0 < a < -0.5 the time lags reach 5 to 20 x r su b/c. It is there- 
fore possible to exploit IR light curves to determine the radial 
brightness distributions and, in turn, the radial dust distribution 
as an alternative to interferometry. We note that the actual shape 
of the CCFs also depend on the auto correlation function of the 
input signal which, however, does not change the overall effects 
and conclusions. 

3.3. Complications when dealing with AGN light curves 

One of the complications when trying to model optical and 
IR light curves of AGN is a possible dependence of variabil- 
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Fig. 4: 2.2//m-8.5/im model cross correlation function (CCF). 
From dark-blue to light-blue the distribution changes from very 
extended (a = 0.0) to very compact (a = -2.0). 



ity on frequency. While the model uses a frequency-integrated 
versi on of dL/L, observations are usually in one specific wave- 
band. Meusing er et al.l d201 ll) report Afa oc A~ 2 , converting to 
A/ v oc v° (i.e. independent of frequency), so that A/ v // v (and 
correspondingly dL v /L v ) depends on the slope of f v . In the near- 
IR the slopes are assumed f v oc y^ 3 and probably redder towards 
the optical and UV (e.g. IZheng et alj 1 1 997t fVanden Berk et alj 
l200T[|Kishiinoto et alJ l2008). This means that A/ v // v in the UV 
is larger than in the optical, and because more energy is radi- 
ated in the UV, any light curve observed in the optical poten- 
tially underestimates the variability amplitude of the absorbed 
emission^ On the other hand, in the process of radiative transfer 
the whole accretion disk spectrum is integrated over all frequen- 
cies, convolved with the dust absorption efficiency, a process that 
probably smoothens out the extreme variability amplitudes. We 
can introduce a constant factor w x = (dL/L)/(dL x /L x ), the vari- 
ability efficiency factor at waveband x, as a free parameter to 
compensate for these complications without requiring further as- 
sumptions on the wavelength-dependence of the variability. In 
this way w x will also include any effects that potentially alter the 
relative variability (e.g. non-variable light from the host), but are 
not included in this simple model. 

An additional complication is a possible change in dust com- 
position with distance from the AGN. Several studies found ev- 
idence that the inner torus is dominated by large graphite grains 
while the outer torus contains a mix of silicates and graphites. 
This leads to comparably small sublimation radii and a change 
of emissivity from the near-IR to mid-I R resulting in a "bump" 
at around 3 urn in many type 1 SEDs (e.g . [Kishimoto et al. 2007; 
iMor et al.ll2009HKishimoto et alj|201 lbh . Therefore, when com- 
paring observed light curves in the near- and mid-IR it may be- 



3 We note that this argument can also be applied to reverberation 
mapping of optical broad Balmer emission lines where the variability 
of ionizing photons is probably larger than the variability observed at 
optical wavelengths. 



come necessary to adjust the absorption efficiencies in Eqs.[3]& 
[6] accordingly. 



4. Modeling the near-IR light curve of NGC 4151 

A comparison between near-IR and mid-IR light curves seem to 
be a promising tool to constrain the brightness profiles of AGN 
dust tori. As of yet mid-IR light curves are very sparse. However, 
as described in Sect. I3.2.U near-IR light curves have different 
peak intensities and tails for compact or extended distributions 
with respect to the AGN variability signal. We may, thus, hope 
that comparing light curves in the optical and near-IR can also 
help decide if the brightness distribution in the torus is compact 
or extended, at least in the inner part of the torus. 

One of the best-monitored AGN in the IR is the b right 
nearby Seyfert 1 galaxy NGC 4151. iKoshida et al l J2009) pre- 
sented optical and near-IR light curves with a very high tempo- 
ral coverage, and the data presented below were extracted from 
their wo rk. They were tak en in the course of the MAGNUM 
project dYoshi et all l2003h over about 2000 days from 2001 
to 2006. We use the reduced and host- and accretion- disk- 
sub tract ed data as shown in F ig. 1 of IKoshida et al.l d2009t) . As 
noted in lMinezaki et al.l (120041) the typical photometric errors for 
NGC 4151 d ata are 0.01 mag (~1 % in flux; but see below for 
more details). Koshida et al. (2009) report that the K-band emis- 
sion at 2.2/im varies similarly to the V-band emission with a 
notable delay that reflects the light-travel time from the accre- 
tion disk to the torus. Near-IR interferometry of the same object 
showed that the /^-band emission is probably coming from a 
confined region at the dust subli mation radius [ Kishimoto et al.l 
(i.e. consistent with a thin ring; 2009b, 2 01 lab . Therefore it 
seems that the object is well-suited to apply our model. 

The AGN variability model has essentially two parameters: 
the brightness distribution power law index a (see Sect. 13.21 ) and 
the variability efficiency factor wy (see Sect. 13.31 ). For the simu- 
lations we first determined the V-to-K time lag by Monte Carlo 
simulations of light curves based on t he observed V- and K - 
fluxes, following the method outlined in Suganu ma et al.l (12006). 
For that we first calculated the structure functions s(tj - tj) = 

(Li<j[f(td - f(tj)fj /N(i < j) of the V- and K-band light 
curves, where /(£,-) is the flux at epoch f, and N(i < j) is the num- 
ber of (/, j) pairs. From s we can determine the typical change 
of flux for a given time interval between observations. Based on 
s we perform Monte-Carlo simulations to interpolate the gaps in 
the observed light curves as follows: First we define a regular 
f-sampling (T\, T2, ■■■) that is finer than the observed one. Then 
linearly interpolated fluxes for this f array are calculated based 
on the observed fluxes. We pick a random epoch T, from the pre- 
defined f-sampling, calculate the distance Tj - f D b s to the near- 
est observed epoch f b s , and determine the corresponding flux 
change A/,(r, - f D b s ) from the structure function s. We then cal- 
culate a random flux change Af from the linearly interpolated 
flux at Tj based on a Gaussian distribution with standard devia- 
tion Af s . The corresponding epoch and flux is added to the array 
of observed fluxes and the procedure repeated until fluxes are 
calculated for all T, in the f-sampling. For a proper comparison 
we use the same f-sampling for the V- and /T-band light curves. 

The Monte Carlo simulations are repeated 50 times to obtain 
an idea of the uncertainties in the light curves introduced by the 
finite sampling of the observations. These simulated light curves 
can then be used to analyze the time lag between the V- and K- 
band emission. We find an average time lag measured over the 
full light curve, excluding the strongest (and widest) peak, of 
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Fig. 5: Observed and interpolated V-band light curve of 
NGC 4151. The blue-filled circles are the observed V data 
points. Error bars are plotted but smaller than the symbols in 
most cases. Based on these observations, 50 interpolated light 
curves have been calculated and shown as gray-dotted lines, re- 
flecting the uncertainty of the "true" light curve. The red-dashed 
line indicates the mean of all these models which has been used 
as the input signal to our model simulations. 



43.8 + 8.5 days, which corresponds to a (sublimation) radius of 
r su b = 0.037 + 0.007 pc, in ag reement with near-IR interferom- 
et al.ll2009btlPott et al.ll2010t [Kishimoto etalj 



etry (Kishimoto 
1201 lah . When using only individual features of the light curves 
the error is generally larger but the resulting time lag is consis- 
tent with the average one. We will discuss consequences of this 
seemingly constant time lag below. 

For the modeling we used the average time lag as an offset 
between the V and K light curves for the model simulations. In a 
next step the simulated V-band light curves have bee used to cal- 
culate a mean V-band light curve. The observed and simulated 
V-band light curves, as well as the mean light curve, is shown 
in Fig. [5] This mean V-band light curve was used as the input 
signal for our torus variability model. 

In Figs.|6]and|7]we compare the observed /T-band light curve 
(blue-filled circles) and the CCFs based on observations (gray- 
dotted lines), respectively, to the best-fitting model. The blue- 
filled circles in Fig. [6] are individua l photometric observations 
extracted from iKoshida et alJ (120091) . Based on the previously 
described Monte Carlo interpolation scheme, we simulated 50 
(nearly continuous) light curves accounting for the gaps in the 
observed light curve and the resulting uncertainty. They are 
shown as gray-dotted lines, as light curves in Fig.|6]and as CCFs 
in Fig. [7] The range of these simulated curves reflect the range of 
curves if the temporal coverage would be infinite. The simulated 
CCFs based on the observations have been normalized so that the 
peak value is unity. Overplotted (red-dashed lines) is the best- 
fitting /T-band model in both plots. The formal fitting was done 
using the light-curve only. In order to asses the quality of the fit 
we used two differe nt versions of t he err or estimates of the ob- 
servation. Although Kos hida et al.l (l2009h do not provide typical 



Fig. 6: Observed and modeled Tif-band light curve of NGC 4151. 
The blue-filled circles are the observed data points. Error bars 
are plotted but smaller than the symbols in most cases. Based on 
these observations, 50 interpolated light curves have been calcu- 
lated and shown as gray-dotted lines, reflecting the uncertainty 
of the "true" light curve. The model light curve is overplotted as 
a red-dashed line. 



errors, an initial fraction of the NGC 4151 jf-b and photometry 
using the same instrumentation was published in Min ezaki et al.l 
(2004) and a photometric error of about 0.01 mag (or about 1% 
in flux) was reported. When using this value for the full data set 
we obtain a reduced xl = 4.2. An alternative approach to esti- 
mate the error in the data comes from the structure function that 
we used in the Monte Carlo interpolation of the light curves. The 
structure function of an AGN is supposed to follow a power law 
from short to long intervals f; - tj until a break at long intervals 
(e.g. ISuganuma et"ai]|2006l) . At short time intervals the intrin- 
sic variation should approach for f, - tj —> 0. In the TT-band 
data we see, however, a flattening of the structure function for 
f,- -tj ^ 6 days, i.e. the (squared) variation of the fluxes becomes 
independent of the observed interval. This is more typical for 
uncertainty in measurements than intrinsic variability. When in- 
terpreted in this way the measurement error is 3.5% and the best 
fit has^-J = 1.2. 

Despite the simplicity of our model, it reproduces the over- 
all peaks and dips of the observations and interpolations (gray 
region) quite well. There are, however, some notable deviations, 
reflected in the moderate xl when assuming 1 % photometric er- 
ror. In particular the observed light curve drops stronger after the 
major outburst at 20 r su b/c leading to a lower general A/(f)//o at 
the peak around 33 r su b/c although the shape of this peak is well 
reproduced. This general flux level is only recovered at the last 
peak after 40 r su b/c. 

The nominal best-fit model has a rather compact brightness 
distribution with a = -1.75. However the probability distribu- 
tion in a parameters is very broad ranging from about -0.5 to 
-2.0 (and probably beyond) with very similar^ values (assum- 
ing 1 % photometric error), meaning that this parameter is poorly 
constrained. A combination of /f-band with longer wavelength 
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Fig. 7: V/TT-band cross-correlation functions of NGC 4151. 
Based on the observed V- and 7^-band light curves, 50 interpo- 
lated light curves for both V- and /T-band have been simulated 
and CCFs calculated accordingly, shown as gray-dotted lines. 
These reflect the observations and the uncertainty in the determi- 
nation of the time lag due to the gaps and non-uniform sampling 
of the data. The CCF based on the light-curve model is overplot- 
ted as a red-dashed line. All CCFs have been normalized to their 
peak value. 



light curves as discussed in Sect. l3.2l is needed to better constrain 
the dust distribution of the hot-dust/sublimation zone emission. 
The efficiency factor wy is much better constrained. In Fig. [8] 
we show the probability distribution for a and wy for both 1 % 
and 3.5% photometric error. In the range of acceptable a-values 
around the best fit we find wy ~ 0.4+q 2 , meaning that only about 
half of the power in the V-band variability is converted into 
variability in the /f-band. Based on the discussion in Sect. 13.31 
one may have expected that the V-band probably underestimates 
the amplitude of the (relative) variability of the integrated AGN 
emission and, therefore, leading to wy > 1 . However a number 
of other factors can play a role in the value of wy, e.g. uncer- 
tainties in host and accretion disk subtraction in either of the 
wavebands, or additional light close to the nucleus that does not 
participate in the variability, at least not on the covered time 
scales (see below for such a scenario). A more extended study 
using a sample of objects and/or additional UV/optical wave- 
bands could help to better understand energy conversion from 
the accretion disk emission into the dust emission and its obser- 
vational caveats. 

Interestingly, the observed change in flux in the V-band 
does not seem to have show n then naively expected change of 
the time lag . Koshidaet al. (2009) already note that a poten- 
tial change in time lag does not scale with (AL) 1 ^ 2 . With a V- 
band flux change of about a factor of 20 from peak to val- 
lejQ, we could have naively expected that the observed time 
lag/sublimation radius changed by a factor of 20'^ 2 = 4.5 (see 



Fig. 8: Probability distribution of our model grid for the radial 
brightness distribution power law index a and the V-band ef- 
ficiency factor wy. The gray-shaded regions with red-dashed 
boundaries show the probabili ty distribution assumin g an error 
of the observations of 1 % from Minezaki et al ] (12004 leading to 
xl - 4.2. The blue-dotted lines assume an error of 3.5% based 
on the analysis of the structure function with axl - 1-2- 



4 This factor seems t o be a little smaller in the light curve of 
IShapovalova et "alli l2008h covering partly the same period. 



Eq.Q). Koshidaet al. (2009) conclude that the time lag changes 
from about 60-70 days in the first 2/3 of the light curve to about 
35-50days in the last 1/3, although at different scaling than 
(AL) 1 / 2 . If such a change were significant, we should have no- 
ticed a shift in the peaks when comparing observed and modeled 
light curves because our models assumes a constant time lag. 
Over about 600days (~ 13.7 r su b/c) the shift between observed 
and modeled light curve would be ~ 6 r su b/c. However, such a 
shift is not seen in Fig. [6] All the peaks and valleys in modeled 
and observed light curve overlap within less than 1 r su b/c. This 
would have also affected the comparison of observed and mod- 
eled CCF, but both are consistent within errors (see Fig.[7]l. 

In order to test even small effects of possible dust destruction 
and reformation, we removed the assumption of a constant time 
lag/sublimation radius and account for dust sublimation once it 
heats over the sublimation temperature (set as a free parameter) 
and instant or delayed reformation of the dust after cooling. This 
changes the sublimation radius and time lag as the variability 
progresses through the torus. However we find that the smallest 
xl is always found for a model with constant sublimation ra- 
dius/time lag, i.e. disfavoring significant dust destruction or ref- 
ormation over the observed time span. This implies that the dust 
can either strongly overheat before it is destroyed, or that the dust 
is efficiently self-shielded. The survival of a dust grain depends 
on the balance between gas pressure and vapor pressure. If par- 
tial gas pressure dominates, then a dust grain is stable; otherwise 
it evaporates. The vapor pressure of dust, p vap oc exp(-l/r), is a 
strong function of the temperature T, while the partial pressure, 
Pgas x T, depends only linearly on T. Therefore, we would ex- 
pect that the lifetime of individual dust grains is short once they 
are heated over the sublimation temperature. As a consequence 
the observed behavior would favor shielding in a locally-dense 
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environment instead of overheating of dust grains. This is consis- 
tent with the idea of a clumpy torus where the dust is confined in 
optically-thick clouds. A change in luminosity first acts on indi- 
vidual clouds which may sublimate part of their dust content but 
can resist longer overall at the same location than smoothly dis- 
tributed dust that is not shielded locally. This scenario may also 
explain the low Wy value: If individual clouds close sublimation 
radius are heated up to the sublimation temperature (and above 
in corresponding equilibrium temperature) and their dust gets 
only gradually sublimated from the surface (e.g. as a "melting 
snowball"), their actual peak temperature would remain essen- 
tially constant, leading to a more or less constant K-hwA flux 
over some time. Hence only cooler clouds would contribute to 
the variability on the same time scales as the incident radiation. 



5. Summary and Conclusions 

We present model simulations of time-variable infrared emis- 
sion from dust as a consequence of variability of the incident 
radiation. For that we first introduce a generalized treatment for 
temperature variations, which can be used for all kind of dusty 
environments. We apply this scheme to a simplified model of 
a (clumpy) dusty torus around AGN and investigate how vari- 
ability of the accretion disk radiation influences the torus emis- 
sion in the near- and mid-IR. The main parameter of this model 
is the radial brightness distribution of the torus that has pre- 
viously been shown to be connected to the radial distribution 
of the dust in the torus. We showed that any variability signal 
in the optical is smoothened stronger if the brightness distribu- 
tion is very extended. While this effect is true for both the near- 
and mid-IR, longer wavelengths show much wider transfer func- 
tions than short wavelengths. The time lags between the optical 
and near-IR emission is mostly representing the light travel time 
from accretion disk to the sublimation radius independent of the 
brightness distribution. For mid-IR wavelengths, however, time 
lags can become very long, up to 10s of r su b/c. The effect that 
the brightness distribution influences the time lags seems to be 
much stronger than any similar effect from inclination (at least 
in type 1 AGN) or d etails of the shape of the inner torus, as 
recently presented by Kawaguchi & Moril d201 ll) . This change 
of lag time from near-IR to mid-IR can be used to quantify the 
brightness distribution of the torus, either by comparing optical 
light curves to near-IR and mid-IR light curves, or by directly 
comparing near-IR to mid-IR light curves. 

Finally, the model has been applied to the optical and near- 
IR light curves of the nearby Seyfert 1 galaxy NGC 4151. We 
show that the simple model can reproduce the overall observed 
variability signal with some deviations in the details. For an even 
better match it may be necessary to use the variability scheme in 
the framework of mo re complex torus model, e. g. by incorpo- 
rating it into CAT3D dHonig & Kishimoto1l2010l) . Nevertheless 
we conclude from our modeling that a single-wavelength near- 
IR light curve is probably not enough to constrain the brightness 
distribution in the torus, requiring at least one other light curve 
at longer wavelengths as discussed before. We found, however, 
that about 40% of the energy in the variability signal in the V- 
band has been converted into ^-band variability. This may be 
explained by a scenario where individual clouds close to the sub- 
limation radius gradually sublimate their dust from the surface 
inward ("melting snowball"), essentially keeping their temper- 
ature constant, and, therefore, do not contribute significantly to 
the variability at the measured time scales. We also note that our 
modeling does not support a change of time lag/sublimation ra- 



dius over the observed light curve epoch in spite of a significant 
change in V-band emission. 
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